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Abstract 

The regenerator is a key component to Stirling cycle machine efficiency. Typical regenerators are 
of sintered fine wires or layers of fine-wire screens. Such porous materials are contained within 
solid-wall casings. Thermal energy exchange between the regenerator and the casing is important to 
cycle performance for the matrix and casing would not have the same axial temperature profile in 
an actual machine. Exchange from one to the other may allow shunting of thermal energy, reducing 
cycle efficiency. In this paper, temperature profiles within the near-wall region of the matrix are 
measured and thermal energy transport, termed thermal dispersion, is inferred. The data show how 
the wall affects thermal transport. Transport normal to the mean flow direction is by conduction 
within the solid and fluid and by advective transport within the matrix. In the near-wall region, both 
may be interrupted from their normal in-core pattern. Solid conduction paths are broken and scales 
of advective transport are damped. An equation is presented which describes this change for a wire 
screen mesh. The near-wall layer typically acts as an insulating layer. This should be considered in 
design or analysis. Effective thermal conductivity within the core is uniform. In-core transverse 
thermal effective conductivity values are compared to direct and indirect measurements reported 
elsewhere and to 3D numerical simulation results, computed previously and reported elsewhere. 
The 3-D CFD model is composed of six cylinders in cross flow, staggered in arrangement to match 
the dimensions and porosity of the matrix used in the experiments. The commercial code FLUENT 
is used to obtain the flow and thermal fields. The thermal dispersion and effective thermal 
conductivities for the matrix are computed from the results. 

Nomenclature 


H Distance of temperature gradient, m 
n Unit tensor 

n Shape parameter 

r Radial position (from the centerline), m 

u Local velocity vector, m/sec 

y Position normal to the wall, m 

y m Dimensionless cross-stream 

effective conductivity 
z Dimensionless streamwise position 

0 Dimensionless temperature 

<js Porosity 


Subscripts: 
d Dispersion 

D Darcy 

eff,yy Effective, cross-stream 
f Of the fluid 

s Of the solid 

tor Tortuosity 

w Wall 

Mathematical Symbols 
Fluctuation 

< > Volume or spatial averaged 


A typical Stirling cycle regenerator geometry is a high-porosity medium made of fine wires inside 
of an annular solid-walled casing. During the cycle of a Stirling machine, the casing and the 
adjacent porous medium may not have the same axial temperature profile, which may lead to heat 



transfer between the porous medium and the casing. This heat transfer may degrade the cycle 
thermal efficiency below that of an ideal engine in which the matrix is isolated from the casing. 
Thus, it is imperative that the mechanism for heat transfer from matrix to casing is well understood. 

In the case of convective heat transfer within a porous medium interrupted by an impermeable wall, 
one can first consider the region that is too far from the wall to be influenced by it. Within this 
region, there is thermal transport due to both 1) hydrodynamic mixing of the fluid within the pores, 
including that by eddies created as a result of flow separating from the solid elements, termed 
thermal dispersion, and 2) conduction through the solid and fluid phases of the porous material. 
Within this region, one can volume average the energy equation over the two phases (solid and 
fluid) of the porous medium. 

The energy equation for the fluid phase is: 
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while the energy equation for the solid phase is: 
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One now can volume average over a Representative Elementary Volume (REV) and, if one assumes 
local thermal equilibrium, a single governing energy equation results, as described by Whitaker 
(1999) and Nakayama and Kuwahara (2000). 


M>c, ),+(!- <fXpc). It) + (pc, ), V • 4(a)(T) = 

V ! (*k f + (1 - 40k. )(T) + V ~ |Mr - k .T> - V '(pc, ), (fa) 


(3) 


The second to the last term in equation (3), the left hand side of equation (4), represents the 
volumetric heat conduction through the phases of the porous medium. It is strongly dependent on 
the structure of the porous medium. This term accounts for the tortuous path that energy must 
follow as it passes through the two phases. It is condensed into another expression by using the 

tortuosity conductivity, k tor ■ 

V-l|(k f T f -k s T s )ds = k t orV{T) (4) 


The first two terms on the right hand side of equation (3), combined with equation (4), are cast as an 
effective stagnant conductivity, k o . 

ko = (<)>k f + (l - <J))k s )n + ktor (5) 


The last term of equation (3) takes the thermal dispersion due to fluid advection into account. It is 
cast in terms of a dispersion conductivity, k d . 



Equations (5) and (6) are often combined to create a total effective conductivity, k e ff . Using it 
results in the following equation: 
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Let us reconsider the porous medium. Within this porous medium there are solid wires and a 
moving fluid. Far from the wall, this fluid will exhibit transport of energy norma! to the streamwise 
direction due to interstitial, or pore scale, mixing. However, as we are volume averaging over the 
fine detail, this value should remain uniform within the core of the porous medium, and will be 

reflected by a single value of the total effective conductivity, keff within the porous medium. 

Consider a region of fluid close to the impermeable boundary, or casing. When one moves radially 
inward from the wall through this region, one passes from a 0% porosity wall to the porosity of the 
core, 90% in the present study (which is typical of a Stirling cycle engine regenerator). Near the 
wall, the uniformity of the matrix and its properties will be disturbed. Even though the edges of the 
screen layers, which are stacked to form the porous matrix, are snugly against the casing wall the 
conduction path is broken (the matrix and the wall are not fused together). This will result in a 
contact resistance at the wall, not present in the core of the matrix. The fluid flow around wires in 
the near vicinity of the wall may also be modified. The flow resistance may be altered as a result of 
the geometry change and its effects on the advective transport of momentum within the near-wall 
porous material. This is probably predominantly due to dampening of eddies near the casing wall. 
We would expect this to be present at distances from the wall that are of the scale of the largest 
eddies within the matrix, which will be a few hydraulic diameters of the matrix. The change in eddy 
transport near the wall will also affect the cross-stream transport of thermal energy over the same 
layer. Also, the contact resistance is high compared to the resistance due to conduction within the 
matrix over the same distance. This is expected when the conductivity of the solid phase is high 
relative to the conductivity of the fluid phase, as is the case with air and stainless steel in the present 
experiment and helium and stainless steel in a representative Stirling engine. Due to these near-wall 
features, one expects the cross stream component of effective conductivity to vary radially from a 
value approaching the molecular conductivity of the fluid very near the wall to the sum of the 
effective molecular conductivity plus dispersion of the in-core matrix. 

The streamwise developments of near-wall momentum and thermal fields resemble boundary layer 
behavior, but with a particular distribution of effective conductivity. Consider a boundary layer 
flow uninfluenced by gravity while assuming slug flow far from the casing (due to the effective 
mixing within the porous medium). As given by Burmeister (1993), the energy equation becomes: 
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If one applies that equation to a steady state, cylindrical, two-dimensional geometry with 
unidirectional flow and a variable cross-stream effective conductivity, k, a negligible pressure work 
term effect, pV • u = 0, negligibly axial conduction, no viscous dissipation, <l> = 0 , and then volume 
average the result, one obtains equation (9), where all length variables are nondimensionalized on 
the cylinder’s radius, r 0 . Measurements indicate that the velocity is nearly uniform, even near the 
wall, in the case of this matrix and the flow conditions of the present experiment (McFadden, 2006). 
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and Ti is the inlet temperature to the test section and z is the dimensionless stream wise position. In 
the following, z = 0 will be the location of onset of heating at the inlet to the test section. 

Consider a ratio of the local cross-stream effective conductivity to the fluid conductivity. 


y m 



(ID 


Recognizing that y m must vary from a value of near unity at the casing wall to some uniform value 
away from the wall and that equation (9) requires that y m be once differentiable, an exponential 
decay form is chosen. This is consistent with the exponentially decaying model proposed by 
Metzger et al. (2004). 



( 12 ) 
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Equations (13) and (14) are the constants which meet our boundary conditions and n will be chosen 
to fit the profile of the radial temperature distribution. 
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This results in a radial relative conductivity distribution of the form shown in Figure 1. Of note is 
that for an n value of 10, the wall-affected region is more than 5 hydraulic diameters thick. Larger 
values of n lead to thinner wall-affected regions. It remains to determine which n value is suitable. 

As mentioned, measurements of the exit plane velocity distribution show for this case that the 
average velocity is uniform across the entire matrix, even in the near-wall region. Equation (6) 
shows that the dispersion conductivity (a component of total effective conductivity) takes into 
account any thermal transport by advection, including velocity fluctuations due to eddies within the 
matrix. This term, when taken with equation (5), results in an effective conductivity that includes 
the thermophysical properties of the porous medium, the geometry of the porous medium, and the 
cross-stream advective transport of thermal energy. 

Equation (9) is solved numerically by using Maple 8® to compute a temperature profile. Some 
example solutions are shown in Figure 2. This model allows for steep drops in temperature near the 
wall with relatively shallow curvature further into the core. For comparison, a profile with a 
uniform k is shown. The boundary conditions applied for this solution are 1) symmetry at the axial 
centerline, 2) uniform inflow temperature, and 3) measured wall axial temperature distribution. 



Figure 2: Temperature profiles computed with two n values, and one with uniform k. 




Previous attempts have been made by Niu et al. (2006) and Niu et al. (2005) to describe the cross 
stream transport terms by direct means, indirect means, and computation, all in the same matrix 
geometry and with the same materials and all under the same flow conditions. Niu et al. (2006) 
described the transport terms by using simultaneous measurements of the three components of 
fluctuating velocity of the flow departing the porous medium. Niu et al. (2005) indirectly 
determined the cross stream thermal dispersion within a porous medium using data taken within a 
spreading thermal wake embedded in the porous medium. Temperature measurements were taken 
inside the porous medium in the same manner as utilized in the present study. Velocity was unform. 


The n value is determined by comparison of computed temperature profiles with measured profiles 
in the matrix. The experiment is in a steady flow that is uniform in temperature and velocity upon 
entry to the matrix. The stainless steel wall of the experiment is heated. Temperature profiles are 
measured at various streamwise locations from the onset of heating to the matrix exit. The wall is 
instrumented to provide the wall thermal boundary condition. Temperature traverses are taken 
radially at discrete streamwise locations. The inflow temperature of the fluid is used as a reference 
point to normalize any inflow temperature fluctuations. Measurements indicate that one may 
assume a uniform velocity throughout the matrix for this case. 

The temperature traverses are made using a 0.05-0.1 mm diameter thermocouple probe that is 
inserted between layers of the porous medium. This probe is inserted until it makes contact with the 
far wall of the test section, as verified by an electrical continuity check between the casing wall and 
the thermocouple bead. The last centimeter of this travel is done at a rate of about 0.1 mm/s so as to 
not bend or damage the thermocouple bead. Forty or more temperature readings are taken at each 
radial location. The traverse starts at the wall and moves toward the core initially by 0.01 mm steps 
using a traversing mechanism. Relative location within the profile is found with 0.001 mm 
precision. As the probe leaves the region near the wall, the step size is increased to reduce the 
quantity of data taken. The total traverse path is 100 mm long, which is past the centerline of the test 
section (95 mm). Figure 3 shows a diagram of the entire test section setup. The test section is a 
317mm long by 190 mm diameter stainless steel casing filled with 200 layers of 158 x 158 (per 
meter) square, welded mesh constructed of round wires of 0.81 mm diameter. Each successive layer 
is rotated 45° relative to its neighbors. The test section porosity is 90%. The test section wall is 
heated over the same region as the probe access region shown in Figure 3. 


Variable Flow Quieting Laminar Flow Flow Quieting Test Section 

Speed Fan Chamber Meter Chamber 



Figure 3: Test Section Setup. Matrix is a cylindrical casing filled with 220 layers of 4 x 4 mesh (4 
wires per inch) welded screens each rotated 45 0 from it predecessor while stacking. 
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Figure 4: Experimental profiles of temperature vs radial position from the wall in mm. 
Temperatures are referenced to the inlet flow temperature. 


Measurements are taken at the axial locations referenced in Figure 4, measured from the inlet of the 
test section, which is also the onset of heating. The temperature difference measurement uncertainty 
is 0.05 °C and the position uncertainty is <0.1 mm. These are based on a 95% confidence interval. 
One can notice the steep drop in temperature at the wall followed by a shallowing curve as one 
moves toward the core, as speculated in the previous section. 

Figure 5 shows that with the proper choice of n, an excellent agreement between the analytic 
solution developed previously and the measured thermal profiles can be obtained. 


Dimensionless 

Temperature 



Figure 5: Dimensionless temperature versus radial distance from the centerline. A representative 
comparison between computed and measured temperature profiles. 






Figure 6: Example of an error plot 

One must consider the accuracy with which the solution is obtained. There are two variables that 
will be solved, the core dispersion and the curvature controlling variable, n. This is done by 
computing a large range of possible values for core y m and n and comparing the experimental 
profiles to the analytical profiles. The best fit establishes y m and n. From this, a plot such as Figure 
6 is generated. The white area shown in Figure 6 is a region of minimum error in y m , n coordinates. 
The error is the sum of the squares of the deviations between the predicted results from the analytic 
solution and the points in the experimental thermal profile. For instance, take one experimental 
profile. A matrix of possible solutions is computed from the analysis over a range of y m and n 
values. Then a comparison is made by subtracting solution values from experimental values for 
each point of the profile. A difference is computed over a range of y m and n values. The difference 
at each data point is squared and summed with the others for that single combination of y m and n to 
get a net error value for each set of y m and n. This is plotted, as shown, for all combinations of y m 
and n. Notice that as n becomes large, meaning that the thickness of the low effective conductivity 
region is very small, the value of y m decreases. Also, as the insulating layer at the wall thickens 
(decreasing n) the relative importance of y m increases. 

The results from this technique are listed in Table 1. In-core values of conductivity are plotted in 
Figure 7 (three points at Pe = 102). For comparison, an independent, indirect measurement of 
effective conductivity in the core of the same matrix was previously made with a different 
technique and was reported by Niu et al. (2005) and Simon and McFadden (2006). This technique is 
a steady state version of one documented by Metzger et al. (2004). These in-core, indirectly- 
measured values are also given in Figure 7. The eddy component of the parameter y m (essentially, 
the slope of the effective in-core conductivity vs. Pe line) was measured previously by a direct 
means and was reported by Niu et al. (2006). These results are also shown in Figure 7. As noted, 
this gives only the slope of the line for the in-core effective conductivity. A stagnant conductivity 
(intercept) of 30 was assumed for the purpose of plotting the line based on these direct 
measurements. 




Figure 7: In-core transverse effective thermal conductivity (summation of transverse stagnant 
conductivity and transverse thermal dispersion) versus Peclet number (for k s /kf = 673 “stainless 
steel” and porosity = 0.9). Compared are the present data, another set of indirect measurements 
taken by an independent technique (Niu et al. 2005), direct measurement results by Niu et al. 2006, 
and computed values from Niu et al. 2005. 


Table 1: Shape Factor and Cross Stream Conductivity for Pe : 


Location 

y,n Core 

n 

152.4 mm 

23 

24 

127 mm 

20 

25 

101.6 mm 

15 

25 


Plotted in Figure 7 are conductivity values computed with CFD. These values come from Niu et al. 
2005, but are briefly discussed below. Two different approaches are typically used in numerical 
modeling of heat transfer in porous media: 1) a high-resolution analysis of the matrix 
microstructure, and 2) a continuum and homogeneous analysis based on the spatially-averaged 
Representative Elementary Volume (REV). In this analysis, the first was done by considering a 3D 
periodic wire arrangement in which each wire is welded with three others to form a square layer. 
These layers are staggered with respect to adjacent layers as they are stacked, as shown in Figure 8. 
When a periodic boundary condition is applied, the array represents a repeating section of the 
matrix. This simulation is a replication of the 4 x 4 mesh used in the experiments discussed above, 
except that it does not have the 45° rotation of each successive layer. The FLUENT 6.2.16 package 
was applied to obtain the numerical solution of the flow and temperature fields inside the matrix. 
This CFD model is believed to be the first 3-D model to represent such a geometry, earlier attempts 





were confined to 2-D arrays or arrays with rectangular wires (Nakayama and Kuwahara, 2000). The 
analysis starts with a periodic velocity condition and a linear inflow temperature profile. After 
solving for the flow and thermal fields, outflow velocity and temperature results are taken back into 
the inflow boundary condition and the solver is reloaded, after verifying proper mass flow rate. This 
is repeated until a periodic thermal field is achieved based on equation (15). It usually takes no 
more than 10 iterations to obtain approximately periodic thermal boundary conditions, 
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wall 
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bulk , inlet. 
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wall 


T —T 

outlet wall 
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bulk, outlet, wall 
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where T bulk 



Tiniet and T ou tiet are averaged values. 


After solving for the temperature distribution inside the 3D cell, the y-direction effective thermal 
conductivity due to thermal dispersion and stagnant conductivity are calculated. The flow-normal 

component of the dispersion, k e ff, yy , is computed with the following equation: 
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Where v is the velocity vector, u , projected in the y direction and H represents the distance over 
which the temperature gradient is imposed. 

A parametric study was carried out to investigate the influence of Peclet number on the thermal 
dispersion conductivity in the 3D cell. Parameter ranges were: Pr = 0.707, 0 < Pe < 1500, <|> = 0.9, 


Solid 



Figure 8: 3D Cell representation for CFD 



and, k s /kf= 673 and 



=3273 (stainless steel). One case, without flow, is used to compute 


the stagnant thermal conductivity. In successive cases, the mass flow rate is varied to compute 
thermal dispersion versus Peclet number. The effective lateral thermal conductivity plotted in 
Figure 7 is calculated by adding thermal dispersion and stagnant thermal conductivity. There is 
general agreement with data from the indirect and direct experimental measurements discussed 
above. The computed thermal dispersion component (the slope on Figure 7) from the calculation is 



0.0047 Pe . The Peclet number length scale is the matrix hydraulic diameter. 


Table 1 shows that if there were no effects of the wall, the cross stream thermal transport would be 
essentially uniform. Our data show this to be generally equivalent to values found by another 
indirect measurement method, direct measurement, and computation. The present results also show 
there are lower values of cross-stream thermal transport near the wall. This is expected since the 
conduction path at the wall is not continuous. A model is given for the radial variation of cross- 
stream thermal transport near an impermeable wall. Its development assumes uniform flow velocity 
throughout the matrix. This assumption is supported for this case by our measurements of matrix 
exit velocity distributions. This variability of near-wall transport will significantly affect the heat 
transfer between the porous medium of a regenerator and the casing wall within which it resides. 
Such variability must be considered in an analysis of a Stirling cycle regenerator. 
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